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Numerical simulations of flow around the F-16XL are presented as a contribution to the 
Cranked Arrow Wing Aerodynamic Project International II (CAWAPI-II). The NASA 
Tetrahedral Unstructured Software System (TetrUSS) is used to perform numerical 
simulations. This CFD suite, developed and maintained by NASA Langley Research Center, 
includes an unstructured grid generation program called VGRID, a postprocessor named 
POSTGRID, and the flow solver USM3D. The CRISP CFD package is utilized to provide 
error estimates and grid adaption for verification of USM3D results. A subsonic high angle- 
of-attack case flight condition (FC) 25 is computed and analyzed. Three turbulence models 
are used in the calculations: the one-equation Spalart-Allmaras (SA), the two-equation 
shear stress transport (SST) and the ke turbulence models. Computational results, and 
surface static pressure profiles are presented and compared with flight data. Solution 
verification is performed using formal grid refinement studies, the solution of Error 
Transport Equations, and adaptive mesh refinement. The current study shows that the 
USM3D solver coupled with CRISP CFD can be used in an engineering environment in 
predicting vortex-flow physics on a complex configuration at flight Reynolds numbers. 


I. Nomenclature 


Cp 

surface pressure coefficient 

a 

angle of attack 

M x 

Rn 

U 

free stream Mach number 

Reynolds number based on reference chord 

free stream reference velocity, ft/sec 

8j 

first cell height, in 

Acronyms 

CAWAP 

Cranked Arrow Wing Aerodynamics 
Program 

NATO 

North Atlantic Treaty Organization 

CAWAPI 

ETE 

Cranked Arrow Wing Aerodynamics 
Program, International 
Error Transport Equation 

RTO 

Research and Technology Organization 
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II. Introduction/Background 

T he Cranked Arrow Wing Aerodynamic Project International (CAW API) utilized the F-16XL aircraft 
as part of a basic research project planned in support of the High Speed Research Program (HSRP). A 
review of the project and how it evolved over the years is given by Lamar & Obara [1]. Flight, wind 
tunnel and computational studies were conducted, and various data sets were generated, analyzed, and 
compared to CFD results [2-5]. CAW API objectives were to validate new methodologies, to evaluate a 
number of predictive methods against available flight test data at high Reynolds numbers, and to check 
the technology readiness level (TRL) of Computational Fluid Dynamics (CFD) codes for a military 
aircraft. Military requirements result in a need for a better understanding of the aircraft characteristics 
before full-scale production. For this purpose, new and existing CFD codes have to be validated, and their 
TRL checked and/or increased. CAWAPI members embraced the idea of engaging in a cooperative 
venture. The benefits from validated CFD codes can enhance the analysis of system performance prior to 
flight, as well as tools to aid in the understanding of unexpected flight behavior. CFD results showed 
good comparison with flight data except for flight condition (FC) 25 and FC 70 [3-5]. FC 25 is a subsonic 
high angle-of-attack case at flow conditions of = 0.242, a = 19.84°, and Rn = 32.22 x 10 6 based on 
reference chord length. FC 70 is a transonic low angle-of-attack case where = 0.97, a = 4.37°, Re = 
88.77 x 10 6 . For FC 25, the worst agreement between computed results and flight data was on the 
outboard wing panel where there were strong interactions between the vortices shed from the inboard 
wing panel, actuator pod, air dam, and crank. The reason for discrepancy was attributed to possible flow 
unsteadiness and turbulence model failure to capture vortex breakdown on the inboard wing panel. The 
CAWAPI- 1 recommendation was to generate better grids, with improved density to capture vortical flow, 
to adapt to flow features and to use higher fidelity physical models. 

The present paper is part of the CAWAPI-2. This is an international cooperation with contributions 
from several organizations. The objective is to find out why CFD fails to simulate these two flight 
condition cases. Relevant questions include: 

1) What is the best flow physics based model or “turbulence model” to employ? 

2) Do we need to adapt for vortical flow generated at these conditions? 

3) Do we have a grid converged simulation? 

4) What is the correct geometry for the transonic case FC 70? 

Due to limited time and resources, the present study focused on the subsonic high angle-of-attack 
case, FC 25. The current study addresses the first three questions noted above, while the fourth concern 
will be fully addressed in future work. Other groups [6-12] within CAWAPI-2 will examine these issues 
as well as whether there is any aero-elasticity issue that affects the fourth question. In the present paper, 
the NASA Tetrahedral Unstructured Software System [13] (TetrUSS) flow solver, named USM3D, is 
used to present detailed comparisons of surface pressures for FC 25 condition. The results produced by 
the one-equation Spalart-Allmaras (SA), the two-equation SST, and the ke turbulence models are 
compared with flight data. An error estimation package developed with funding from a NASA Small 
Business Innovation Research (SBIR) project is used to analyze the CFD results and provides guidance to 
the limit of turbulence models used. Grid adaption and uncertainty quantification using the commercial 
code CRISP CFD was conducted and is also presented. 

III. Computational Approach 

In this section, details of the computational grids and the numerical approach used in this study are 
presented. Various convergence criteria adopted to ensure solution convergence are described. 

USM3D Flow Solver 

USM3D is a tetrahedral cell-centered, finite volume, Euler and Navier-Stokes flow solver. Time 
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integration follows the implicit point Gauss-Seidel algorithm, explicit Runge-Kutta approach, and local 
time stepping for convergence acceleration. The solver scheme allows various options for computing the 
in viscid flux quantities across each cell face. These include Roe’s flux-difference splitting (FDS), 
advection upstream splitting method (AUSM), flux vector splitting (FVS), and the Harten, Lax, and van 
Leer with contact restoration (HLLC) scheme. MinMod and Superbee flux limiters have been 
incorporated in the flow solver to smooth out the flow discontinuities due to shock waves. Turbulence 
models include Spalart-Allmaras (SA) [14] with and without the wall function, shear stress transport 
(SST) [15] and ke.[16]. 

Grid Generation 

The surface triangulations along with the field tetrahedral volume grids were generated using the 
GRIDTOOL and VGRID software developed at NASA Langley Research Center (LaRC). A rectangular 
box that encompasses the vehicle is typically used to define the computational domain far-field 
boundaries. Each face of this rectangular box is located several body lengths away from the configuration 
in the upstream, radial and downstream direction. The farfield boundaries were located at 100 times the 
mean aerodynamic chord of the wing. As a general practice, each final converged solution is analyzed to 
insure that the laminar sub-layer, has been grid resolved and that the average y + is less than 1. Four 
different grids (extra coarse, coarse, medium, and fine) have been generated to facilitate grid-converged 
solutions. Table 1 provides the average normal spacing to viscous wall in terms of y + coordinate, the 
values of the first cell height 6 j5 the total number of cells in the viscous layer, and the total number of 
surface elements for different grids. 

Table 1: Grids and computational resources. 


Grid 

Number of 
cells 

First cell 
height, 5 ls in 

Y+ average 

Number of 
viscous cells 

Number of 
surface triangles 

Number of 
iterations 

Number 
of CPU 

Coarse 

19,370,847 

0.0009225 

1.108 

15,193,683 

256,846 

50,000 

288 

Medium 

62,473,588 

0.000615 

0.757 

40,200,735 

640,850 

50,000 

720 

Fine 

143,034,292 

0.00041 

0.516 

83,173,422 

1,265,352 

50,000 

1200 


Error Quantification 

Typical grid refinement studies, while effective in quantifying the discretization error inherent in any 
CFD solution, are often intensive in terms of time, labor, and CPU cost. As a supplement to such studies, 
solution verification is performed using an Error Transport Equation (ETE) solver contained in the CRISP 
CFD package. The ETE method [17-19] has been developed for use with various unstructured grid 
solvers, including USM3D. The basic premise of the ETE is that discretization errors may be generated 
in one location of the domain and convected elsewhere as erroneous waves propagated by the 
characteristics of the Euler/Navier-Stokes equations, where they become apparent and manifest as 
interpolation errors. The source term of the ETE is formed from the leading terms of the truncation error 
associated with the inviscid flux scheme. Once the ETE solution is obtained, errors in any number of 
local or integrated quantities may be examined [18]. 

The inclusion of numerical error bars as part of post-processing and analysis offers new possibilities 
for how CFD results can be interpreted and used. If the increments provided by grid refinement can be 
reliably predicted, a higher level of confidence may be placed in lower fidelity meshes. This would 
potentially preclude the need for fine grid analysis if the numerical errors predicted for a medium 
resolution mesh were judged to be acceptable. In addition, if grid-induced errors are bounded and 
experimental data lies outside numerical error bars, it would likely indicate a deficiency in physical 
modeling. Such information would be of considerable benefit to the CFD analyst in evaluating 
simulation results. 
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Adaptive Mesh Re finement and Coarsening 

CRISP CFD contains a parallel mesh adaptation and quality improvement capability for three- 
dimensional mixed-element unstructured meshes [20-21]. Meshes comprised of tetrahedral, prismatic, 
and hexahedral regions may be readily modified to generate more accurate flow solutions through local 
refinement and coarsening. Tetrahedral refinement and coarsening are accomplished using Delaunay 
refinement and edge collapse methods, while cell subdivisions are employed for viscous tetrahedra, 
pyramids, prisms, and hexahedral cells. Upon completing the mesh modifications, the solution is 
interpolated onto the adapted mesh using a grid-transparent procedure based on nearby point clouds [22]. 

An estimate of the solution error is used to drive mesh adaptation. For steady state applications, the 
default method currently employed [21] is an interpolation error, based on forming a higher order 
approximation of the solution at each mesh point using a least squares approach. The difference between 
the higher order reconstruction from incident nodes and the current solution forms the error indicator. If 
the current mesh is sufficiently fine to support the spatial variation in the solution, the estimated error will 
be low, allowing coarsening to take place. Conversely, a high degree of error indicates that additional 
refinement is needed. This approach has proven successful in a variety of applications and is capable of 
detecting shear layers, separation, vortical flows, and weak gradients in coarse regions, as well as shocks 
and expansions, with no tuning required of the user . If an Error Transport Equation solution is available, 
it is also possible to use the residual of the ETE to drive mesh refinement [18]. 

Solution Development and Convergence 

All solutions are computed using Roe’s flux difference splitting method. Flow is assumed to be fully 
turbulent. Three turbulence models, the SA, the SST and ke, will be used in the present computations. 
Mean flow and turbulence model equations were solved in a decoupled fashion. In the precursor studies, 
not reported herein, it was established that the solutions were insensitive to the manner in which they 
were advanced in time (steady-state versus 2 nd order time-accurate) and to the prescription of an initial 
state (freestream condition versus converged solution at a lower angle-of-attack). Based on these findings 
all solutions were obtained assuming steady state using the implicit backward-Euler scheme. 

In general, solution convergence is evaluated by monitoring variations in all six longitudinal and 
lateral force and moment aerodynamic coefficients. A solution was considered converged when 
fluctuations in these coefficients were reduced to less than 0.5% of their respective average values 
calculated over the previous 2000 iterations. Additionally, all solutions were required to have the 
cumulative mean flow residual errors reduced by at least three orders of magnitude and turbulence model 
residual errors reduced by at least two orders of magnitude. Each case completed 10,000 first-order 
iterations followed by 40,000 second-order iterations for a total of 50,000 iterations. A sample of 
convergence history is shown in Fig. la. The solution achieved 5-orders of magnitude reduction and all 
six force and moment coefficients were converged. Switching from a first-order to a second-order scheme 
caused the spike in the residual at 10,000 iterations. A NAS supercomputer was used to compute all of the 
cases. A typical case was completed in less than 12 hours. 

IV. Results and Discussions 

The subsonic high angle-of-attack flight condition (FC 25) is simulated for the flow around the F- 
16XL. Computational results, surface static pressures, are presented and compared with flight data. FC 25 
has a free stream Mach number of 0.242, an angle of attack of 19.4 degrees, and Reynolds number equal 
to 32.22 x 10 6 based on reference chord length. Vortical-flow phenomena over the upper surface, is 
expected. The F-16XL has two leading edge sweep angles, 70° for the inboard wing panel and 50° for the 
outboard wing panel. It should be noted that all simulations were conducted as steady RANS calculations 
and this approach might mask unsteadiness effects. 

The F16XL surface pressure distribution is shown in Figure lb. The flow structure consists of 
several vortices: i) an inboard vortex originating from the wing apex, ii) an outboard primary vortex 
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originating from the wing leading edge outboard of the crank, and iii) an airdam vortex. All solutions 
exhibit convergence to 4-orders of magnitude. 

Figure lc shows the surface streamlines. The convergence-divergence of the surface streamlines lines 
shows the location of the flow separation and reattachment. Flight data were extracted at different 
fuselage stations (FS) and butt line (BL) cuts. Surface pressure measurements are extracted at various FS 
and BL locations for comparison with flight test data. The black lines in figure Id shows the seven FS 
locations, 185, 300, 337.5, 375, 407.5, 450 and 492.5. Figure le shows the seven BL locations: 55, 70, 80, 
95, 105, 153.5 and 184.5. BL 184.5 is highlighted in red as it presents a great challenge for turbulence 
model results to match flight data and will be discussed in this section. The primary vortex starts at FS 
185 and the secondary vortex starts at around FS 407.5. There are strong interactions in the vortex system 
over the F-16XL. The inner- wing secondary vortex strikes the airdam followed by an airdam vortex 
emanating downstream. Figures la and lb show the outer-wing primary vortex forming at the crank and 
traveling downstream and outboard towards the tip and the missile store. 

In the following sections, we demonstrate the effects of grid refinement and turbulence models in 
predicting the locations and strengths of these vortices. During the last 10 years of NASA’s involvement 
with CAW API, several grids were produced to understand the aerodynamic characteristics of all the flight 
conditions for the F16XL configuration. Table 1 above lists the coarse, medium and fine grids, that were 
used in the present study. We utilized these grids for turbulence models comparisons, grid comparisons, 
grid adaption and error estimates. The 62 million cell “medium” grid was used as a baseline resolution 
mesh for the assessment of turbulence modeling effects. 

Turbulence Model Comparisons 

Data were extracted at seven BL locations: 55, 70, 80, 95, 105, 153.5 and 184.5 (figure lc). Seven FS 
locations, 185, 300, 337.5, 375, 407.5, 450 and 492.5 (figure Id) were also available for comparison. The 
pressure coefficient data from the well-established one-equation SA turbulence model, the SST model and 
the ks two-equation model were compared, with FC 25 data and presented in Figure 2. In general, all 
three turbulence models compared well with flight data on the lower surface of the wing however, it is 
clear that ks provides the best comparisons with upper surface pressure data as shown Figure 2. Figure 
2g, is the only BL cut that all three turbulence models fails to predict the first 30% of the chord. This is 
the outboard BL cut along the wingspan with the smallest chord length, shown as a red line in figure Id. 
For the other BL and FS cuts, kc results compared well with the flight data. On the other hand, both SA 
and SST gave similar results and over-predict the data similar to the results shown in figure 2f. 

The surface pressure coefficient from ks and SA turbulence models are shown in figure 3. The ks 
results are less diffusive than SA results. The streamlines, at BL 105, superimposed on the total pressure 
contours are shown in figure 4 to envision the discrepancies between the results of SA and kc turbulence 
models. The pressure results from kc slowly increased in the first 20% length of the wing (see figure 2e) 
following the flight data. Also, the kc total pressure is much smaller compared with SA results. This 
causes the streamlines to be pushed down toward the wing surface causing the pressure to slowly increase 
as compared with the streamlines from SA that completely missed the wing surface. 

In the following sections, having performed model validation, solution verification exercises 
consisting of grid comparisons, error-estimates and attempts to adapt the grid are discussed in detail. 

Grid Comparisons 

Four different grids (extra coarse, coarse, medium, and fine) have been used to facilitate grid- 
converged solutions. The extra coarse grid was the original grid used in reference [4]. Three new grids 
were generated for the F-16XL; coarse, medium and fine. The guidelines established by the Drag 
Prediction Workshop [23] and the High Lift Prediction Workshop [24], to address solution sensitivity to 
grid refinement, were adopted in the grid generation. The prior computations for FC 25, reference [4], on 
the extra coarse grid, provided the reference data for the present analysis and, in particular, for generating 
the properly sized coarse, medium, and fine grid levels. Figure 5 shows the representative computational 
grids on the surface and the near-field view for typical cross sectional cuts of the coarse, medium, and 
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fine grid. Table 1 above provides the values of 6 l5 the total number of cells in the viscous layer, and the 
total number of surface elements for different grids. A factor of 2/3 was used to size the value of 8 l5 for 
each consecutive grid refinement. Furthermore, a factor of 1.26 was used to construct the viscous layers 
(Sj) in the radial direction by setting the grid growth rate parameters R1 =0.16 and R2 = 0.04 for all the 
grids used in the present investigation. 

Sj= S 7 (1+Rl(l+R2) 7 _ 1 y -1 

In this equation, the variable 8, is the first cell height next to the solid surface in the normal direction. This 
grid system is suitable for cell-centered flow solvers. 

In an effort to improve the grid quality, a grid optimization technique [23] available in CRISP CFD 
was employed to smooth out any local grid skewness or large variation of the cell volumes before the 
launching of the flow solver. This activity was shown to stabilize and enhance the solution convergence 
characteristics. In addition, as shown in Table 1, there was roughly the expected linear relationship 
between the total number cells (N) and the corresponding computer resource requirements. As a result, 
the medium grid computations required approximately 3.2 times the CPU hours to converge the solution 
as those needed for the coarse grid. This similarity in computer usage remained approximately the same 
between the medium and the fine grid. 

The original grid from CAWAPI I study had a 10.2 million-cells grid with a y+ average of 1.2. The 
grid did not follow the criteria of the Drag Prediction Workshop and is not shown in Table 1. More details 
about this grid can be found in reference [4]. Figure 5 shows the computational grid distribution on the 
surface of the F16XL and two sample cuts at FS 185 and BL 184.5 for the coarse, medium, and fine grids. 
There were 15 components, with different sizes and shapes over the F16XL configuration. Such 
geometrical complexity presented a challenge for the grid generation and the flow solver to properly 
capture the associated flow physics and the resulting aerodynamic properties of the vehicle. 

Figure 6 shows comparison between the computed surface pressures from rom the CAWAPI- 1 grid, 
the coarse, medium, and fine grids, and the measured surface pressures at selected stations BL 70, BL 
105, BL 184.5, and FS 300. The surface pressures at BL 70, BL 105, and FS 300 show little effect of grid 
refinement. At these three stations, coarse, medium and fine grid surface pressures are almost identical 
and grid convergence was achieved. In general all computed results are in good agreement with the 
measured data. This is also evident for all the other stations that are not shown in Figure 6. Grid 
resolution effect can be clearly seen further outboard on the main wing at station BL 184.5. Even fine grid 
surface pressures show poor agreement with the measured surface pressures at this station and grid 
convergence was not achieved. The reason for the discrepancy might be due to ks turbulence model 
failure to predict vortex break down and the onset of separation on the outer wing. Vortex breakdown was 
expected but since steady RANS calculations were performed it might have masked the development of 
unsteady separation. To further investigate this discrepancy, higher fidelity approaches (URANS, LES, or 
DES simulations) must be used. 

Error Estimates 

An example of using error quantification for the F-16XL configuration is discussed in this section. 
The results were used to predict the error bar for surface pressure at different stations along the 
configuration. The solution of the error equations will require computer time equivalent to solving the 
main equations. Once the error equations are solved, it is possible to generate error data for most 
aerodynamic quantities of interest. It is expected that these errors would reduce with mesh refinement as 
the vortex resolution is improved. Furthermore, the errors predicted on a coarse or medium resolution 
mesh should contain the fine grid solution. 

Figure 7 presents predicted errors in pressure coefficient due to grid resolution at selected butt line 
and fuselage station locations for the coarse 19M-cell grid. The medium resolution 62M cell results are 
included for comparison. At nearly each station shown, the predicted errors from the solution of the ETE 
span the observed differences between coarse and medium grid solutions. This provides a quantitative 
assessment of the increment expected from refining the mesh, and serves to verify the computation and 


6 

American Institute of Aeronautics and Astronautics 


AIAA 52 nd Aerospace Sciences Meeting 

National Harbor, MD Invited Paper, Special Session on CAW API-2 Aerodynamics 

bound the error due to discretization. Even the vortex reattachment at FS 185, seen in Figure 7e, is nearly 
accounted for by the predicted errors. 

A similar error analysis is depicted in Figure 8 for the 62M-cell medium resolution mesh. Results 
from the 143M cell fine grid are included for comparison. The observed solution increment is well 
captured by the predicted error bars for each of the four butt line locations considered in Figure 8a-8d. 
The fuselage station comparisons indicate areas where the predicted errors fall short of properly 
containing the fine grid results, as evident in Figure 8e, 8g, and 8h, where the inner and outer primary 
vortices and the dam vortex are better resolved. One explanation for the discrepancy is that the ETE is 
solved on the same discretization as the flow solution itself, and thus is subject to similar numerical 
diffusion issues. Thus it is important to match the limiter, gradient evaluation, flux scheme, etc. between 
the flow solver and the ETE solver to the greatest extent possible. It should also be noted that while not 
present in this example, the ETE cannot account for flow features present on one grid and absent in 
another. That is, the simulations must produce self-similar solutions to place confidence in the solution 
verification provided by the error analysis. 

Figure 9 depicts simultaneous visualization of pressure coefficient and its associated error using the 
“error bubble” technique [18], in which spherical glyphs are colored by the selected variable and sized 
according to the magnitude of the error. In the limit, the glyphs vanish in regions of negligible error. 
This technique allows grid-induced errors to be readily rendered and assessed for complex flows. In 
Figure 9 it is apparent that much of the error due to grid resolution is associated with the leading edge 
vortices. Local grid adaptation holds the potential to reduce and minimize such errors and will be 
explored in the following section. 

Adaptive Grid 

In the current study, two approaches were used in the grid adaption. The first was the default solution 
error estimation method implemented in CRISP CFD, which is based on forming a higher order 
approximation of the solution at each mesh point using a least squares approach [21]. If the current mesh 
was sufficiently fine to support the spatial variation in the solution, the estimated error was low, allowing 
coarsening to take place. Conversely, a high degree of error indicated additional refinement was needed. 
The second approach used information from the ETE solver to drive the mesh adaptation process [18]. 
This method of adaptation would not necessarily refine the flow features. It should be noted that since the 
criterion for coarsening was the same as that used for refinement, keying adaptation to a single variable 
might have removed cells in regions where other variables may have exhibited significant variation. One 
also must be cautious with any coarsening when using error sources. This is an area that still needs to be 
explored. 

CRISP CFD was used to evaluate the flow field variation, refine, coarsen, and optimize the mesh in a 
localized manner, and also produced a new, adapted grid and interpolated solution file for the selected 
flow solver. USM3D was then restarted on the newly adapted grid and ran until fully converged. This 
constituted one cycle of adaption. Ideally, mesh adaptation allows one to obtain more accurate, grid- 
converged solutions in an efficient manner. Given a converged solution on an initial mesh, subsequent 
solutions on locally adapted meshes do not require as much computational expense to re-establish a 
converged state. In steady state simulations, adaptation may be repeated as often as desired, until grid 
convergence is achieved or solution errors are below a defined tolerance. 

In this study, three adaption cycles were conducted. Table 2 shows grid size after each adaption cycle. 
Grid sizes increased from 19 million cells on the original coarse grid to 24 million cells for the final 
adapted mesh. The degree of refinement was not aggressive, and in the interest of retaining practical mesh 
sizes, the viscous cell layers were not modified. Viscous cell subdivision [20] tends to increase the 
number of cells dramatically, as cells are split into 2, 4, or 8 cells according to predefined patterns. As a 
consequence of freezing the viscous layers, the surface mesh was not modified in the present calculations. 
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Table 2: Grids used in the mesh adaptation studies. 


Grid 

Number of cells 

Original Mesh 

19,370,847 

1 st Adaptation 

20,189,477 

2 nd Adaptation 

21,553,095 

3 rd Adaptation 

24,299,853 

Adaptation Using ETE Sources 

25,420,230 


Figure 10 presents comparison between Cp distribution on the different adapted grid levels and the 
solution on the coarse and medium grids at selected FS and BL locations. Figure 10 shows slight 
variations between all three adapted grid levels and the 19M cell coarse grid. As the mesh was adaptively 
refined, it is evident that the solution is moving towards a grid converged result. Figure 10a shows that 
the adapted solution did not capture Cp at the wing leading edge. However adapted solution captured Cp 
peak at x/c of 0.19. The adaption outperformed the medium grid in capturing expansion over the wing 
between x/c of 0.19 and 0.3. The adapted and the medium grid solutions are in close agreement for x/c 
>0.3 yet higher than flight data for x/c > 0.6. Figure 10b shows adapted solution is closer to flight data 
than the medium grid solution, yet still failed to capture the Cp distribution on the leading part of the 
wing. As the mesh was adaptively refined, the peak around 2y/b of 0.8, at station FS 185, moved toward 
medium grid solution as shown in figure 10c. Figure lOd shows that the adaptive solutions, at FS 300, 
were grid converged and in close agreement with medium grid results. FS 300 is one of the stations where 
the computed results are grid converged as shown earlier in Figure 6d comparing coarse, medium and fine 
results. The differences observed between the medium grid results and the adapted mesh solutions may be 
attributed to the fact that the surface mesh was not modified in the adaptive grid sequence; hence 
resolution of the vortex reattachment was not enhanced. It is anticipated that allowing the viscous cell 
layers to be subdivided would alleviate this issue, and this may be explored in future work. 

V. Conclusion 

The results for the F-16XL configuration obtained using the TetrUSS flow analysis and CRISP CFD 
systems are presented. The flow analysis was conducted for a subsonic, high angle-of-attack case, where 
Moo = 0.242, a = 19.84°, Re = 32.22 x 10 6 based on the reference chord length. The configuration was the 
focus of the Cranked Arrow Wing Aerodynamic Project International II (CAWAPI-2). All the computed 
results were based on the assumption that the flow was fully turbulent and steady over the entire vehicle 
surfaces. Turbulence model sensitivity of the numerical results was examined using the SA, SST and ks 
turbulence models. Computed surface pressures were compared with corresponding flight data. In 
general, the ks results compared well with the flight data, while both SA and SST gave similar results that 
over-predicted the data. It was observed that the accuracy of the computed surface pressure, from SA, 
SST and ks turbulence models, was degraded closer to the wing tip region. Surface pressures at BL 184.5 
span were substantially under-predicted over most of the main wing cut. It is surmised either a higher 
fidelity turbulence model or detached eddy simulation (DES) may improve the accuracy in resolving the 
wing tip vortex and localized flow separation on the outer wing. Grid convergence study was conducted 
for ks turbulence model and it showed that the solution was grid converged. 

Solution verification was approached using two methods, prediction of the discretization errors 
using the solution of Error Transport Equations (ETE), and local adaptive mesh refinement. The ETE 
solutions were useful in quantifying the error due to grid resolution on coarse and medium level meshes. 
Including the numerical error bars allow for the discrepancies with the test data to be better understood as 
due to physical modeling issues vs. resolution issues. Adaptive refinement reduced the resolution error 
and confirmed that a grid converged result was achieved. Adaptive refinement on the coarse grid captured 
medium grid results. The areas for future investigation include viscous layers and surface mesh 
refinement, the use of error prediction and adaptive refinement for the transonic FC 70 case, and gaining 
additional experience with the ETE source criterion for adapting the grid. 


8 

American Institute of Aeronautics and Astronautics 


AIAA 52 nd Aerospace Sciences Meeting 
National Harbor, MD 


Invited Paper, Special Session on CAW API-2 Aerodynamics 


Acknowledgments 

The authors would like to thank Susan Cliff and Michael Aftosmis from NASA Ames for the fruitful 
discussions during the course of this work. 


VI. References 


[1] Obara, C. J., and Lamar, J. E., “Overview of the Cranked-Arrow Wing Aerodynamics Project International,” Journal of 
Aircraft, Yol. 46, No. 2, 2009, pp. 355-368. 

[2] Lamar, J. E., Obara, C. J., Fisher, B. D., and Fisher, D. F., “Flight, Wind-Tunnel, and Computational Fluid Dynamics 
Comparison for Cranked Arrow Wing (F-16XL-1) at Subsonic and Transonic Speeds,” NASA TP 2001-210629, Feb. 2001. 

[3] Boelens, O. J., Badcock, K. J., Elmilgui, A., Abdol-Hamid, K. S., and Massey, S. J., “Comparison of Measured and Block 
Structured Simulations for the F-16XL Aircraft,” Journal of Aircraft, Yol. 46, No. 2, 2009, pp. 377-384. 

[4] Gortz, S., Jirasek, A., Morton, S. A., McDaniel, D. R., Cummings, R. M., Lamar, J. E., and Abdol-Hamid, K. S., “Standard 
Unstructured Grid Solutions for CAW API F-16XL,” Journal of Aircraft, Vol. 46, No. 2, 2009, pp. 385-408. 

[5] Rizzi, A., Jirasek, A., Lamar, J. E., Crippa, S., Badcock, K. J., Boelens, O. J., “Lessons Learned from Numerical Simulations 
of the F-16XL Aircraft at Flight Conditions,” Journal of Aircraft, Vol. 46, No. 2, 2009, pp. 423-441. 

[6] James M. Luckring; Arthur W. Rizzi, “Toward Improved CFD Predictions of Slender Airframe Aerodynamics Using the F- 
16XL Aircraft (CAW API-2):, to be presented in 52nd AIAA Aerospace Sciences Meeting, National Harbor, Maryland, January 
13-17,2014. 

[7] Maximilian M. Tomac; Adam Jirasek; Arthur W. Rizzi, “Factors Influencing Accurate Shock Vortex Interaction Prediction 
on F-16XL Aircraft”, to be presented in 52nd AIAA Aerospace Sciences Meeting, National Harbor, Maryland, January 13-17, 
2014. 

[8] Adam Jirasek; Luca Cavagna; Luca Riccobene; Sergio Ricci; Arthur W. Rizzi, “Estimating first-order aeroelastic effects on 
F-16XL using NeoC ASS/Edge analysis”, to be presented in 52nd AIAA Aerospace Sciences Meeting, National Harbor, 
Maryland, January 13-17, 2014. 

[9] Adam Jirasek; Maximilian M. Tomac; Arthur W. Rizzi, “Improved Predictions for high Angle of Attack Flow around the F- 
16XL Aircraft”, to be presented in 52nd AIAA Aerospace Sciences Meeting, National Harbor, Maryland, January 13-17, 2014. 

[10] Okko J. Boelens; James M. Luckring, “Flow Analysis of the F-16XL Aircraft (CAW API-2) At Transonic Flow Conditions”, 
to be presented in 52nd AIAA Aerospace Sciences Meeting, National Harbor, Maryland, January 13-17, 2014. 

[11] Stephan M. Hitzel, “Vortex Flows of the F-16XL Configuration - CAWAPI-II Free-Flight Simulations”, to be presented in 
52nd AIAA Aerospace Sciences Meeting, National Harbor, Maryland, January 13-17, 2014. 

[12] Arthur W. Rizzi ,”What was Learned in Predicting Slender Airframe Aerodynamics with the F16-XL Aircraft”, to be 
presented in 52nd AIAA Aerospace Sciences Meeting, National Harbor, Maryland, January 13-17, 2014. 

[13] TetrUSS Web page: http://tetruss.larc.nasa.gov/usm3d/index.html, May 2013. 

[14] Spalart, P. R., Allmaras, S.R.: A One-Equation Turbulence Model for Aerodynamic Flows. AIAA Paper No. 1992-0439. 

[15] Menter, F. R., “Improved Two-Equation k-omega Turbulence Models for Aerodynamic Flows”, NASA TM 103975, October 
1992. 

[16] Jones, W. P., and Launder, B. E., “The Prediction of Laminarization With a Two-Equation Model of Turbulence,” Int. J. 
Heat & Mass Transf., vol. 15, no. 2, Feb. 1972, pp. 301-314. 

[17] Cavallo, P. A., and Sinha, N., “Error Quantification for Computational Aerodynamics Using an Error Transport Equation,” 
Journal of Aircraft, Vol. 44, No. 6, 2007, pp. 1954-1963. 

[18] Cavallo, P. A., O’Gara, M. R., and Schikore, D. R., “Software System for Prediction, Visualization, Analysis, and Reduction 
of Errors In CFD Simulations,” AIAA Paper 2009-3649, 19th AIAA Computational Fluid Dynamics Conference, San Antonio, 
TX, June 22-25, 2009. 

[19] Cavallo, P. A., O’Gara, M. R., Feldman, G. M., and Liu, Z., “Unified Error Transport Equation Solver for Solution 
Verification on Unstructured Grids,” AIAA Paper 2012-3345, 42nd Fluid Dynamics Conference, New Orleans, LA, June 25-28, 
2012. 

[20] Cavallo, P. A., and Grismer, M. J., “A Parallel Adaptation Package for Three-Dimensional Mixed-Element Unstructured 
Meshes,” Journal of Aerospace Computing, Information, and Communication, Vol. 2, No. 11, 2005, pp. 433-451. 

[21] Cavallo, P. A., Sinha, N., and Feldman, G. M., “Parallel Mesh Adaptation Method for Moving Body Applications,” AIAA 
Journal, Vol. 43, No. 9, 2005, pp. 1937-1945. 

[22] Baker, T. J., “Interpolation from a Cloud of Points,” Proceedings of the 12 th International Meshing Roundtable , Santa Fe, 
NM, September 14-17, 2003. 

[23] Baker, T. J., and Vassberg, J. C., "Tetrahedral Mesh Generation and Optimization," Proc. 6 th International Conference on 
Numerical Grid Generation, ISGG, pp. 337-349, 1998. 

[24] 5th AIAA CFD Drag Prediction Workshop, AIAA, New Orleans, LA, http://aaac.larc.nasa.gOv/tsab/cfdlarc/aiaa-dpw// 
[retrieved 21 November 2013]. 


9 

American Institute of Aeronautics and Astronautics 


AIAA 52 nd Aerospace Sciences Meeting 
National Harbor, MD 


Invited Paper, Special Session on CAW API-2 Aerodynamics 


[25] 2nd AIAA CFD High Lift Prediction Workshop (HiLiftPW-2), AIAA, San Diego, CA, http://hiliftpw.larc.nasa.gov/ 
[retrieved 21 November 2013]. 



a) Convergence history 



Inner Primary Vortex 


Air-Dam Vortex 


Outer Primary Vortex 


b) Surface pressure distribution 



c) Surface streamlines 



d) Fuselage stations 



BL 55 
BL 70 
BL 80 

BL 95 
BL 105 


e) Butt lines 


Figure 1. USM3D solution of the F16XL at Moo = 0.242, a = 19.84°, and Rn = 32.2xl0 6 . 
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(a) BL 55 


(b) BL 70 




(c) BL 80 (d) BL 95 




(e) BL 105 (f) BL 153.5 

Figure 2. Comparisons of computed turbulence model results and measured flight C p at FC 25 

(M„ = 0.242, a = 19.84°, and Rn= 32.2 x 10 6 ). 
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(g) BL 184.5 
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(k) FS 375 




Figure 2. Continued. 
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Figure 2. Concluded. 



Figure 3. Comparisons of upper surface C p using SA and ks models at FC 25 

°, and Rn= 32.2 x 10 6 ). 




b) BL 105, ke 


Figure 4. Computed total pressure and streamlines at BL 105 for FC 25 using SA and 
ke models (M® = 0.242, a = 19.84°, and Rn= 32.2 X 10 6 ). 
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g) BL =184.5 Coarse Grid Cut h) BL =184.5 Medium Grid Cut i) BL =184.5 Fine Grid Cut 

Figure 5. Computational grid distribution on the surface and at FS 185 and BL 184.5 

for the coarse, medium, and fine grids. 
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(c) BL 184.5 



Figure 6. Comparisons of flight data with ks turbulence model Cp on different grid 
levels (Moo = 0.242, a = 19.84°, Rn= 32.2 X 10 6 ). 
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(b) BL 70 



(c) BL 80 





Figure 7. Computed error estimates on 19M-cell mesh compared with 62M-cell 
medium mesh results (M* = 0.242, a = 19.84°, Rn= 32.2 X 10 6 , ks turbulence model). 
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(h) FS 450 


Figure 7. concluded. 
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(c) BL 80 (d) BL 153.5 

Figure 8. Computed error estimates on 62M-cell mesh compared with 143M-cell fine 
mesh results (M* = 0.242, a = 19.84°, Rn= 32.2 X 10 6 , ks turbulence model). 
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(h) FS 450 

Figure 8. concluded. 



Figure 9. Bubble plot visualization of errors in pressure coefficient for 19M-cell coarse mesh. 
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(c) FS 185 (d) FS 300 

Figure 10. Comparisons of computed Cp at different grid adapted levels and measured 
flight data (M*, = 0.242, a = 19.84°, Rn= 32.2 X 10 6 , ke turbulence model). 
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